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Abstract 

Previously we have developed a free boundary model for local ther- 
mal coagulation induced by laser light absorption when the tissue region 
affected directly by laser light is sufficiently small and heat diffusion into 
the surrounding tissue governs the necrosis growth. In the present paper 
surveying the obtained results we state the point of view on the necrosis 
formation under these conditions as the basis of an individual laser ther- 
apy mode exhibiting specific properties. In particular, roughly speaking, 
the size of the resulting necrosis domain is determined by the physical 
characteristics of the tissue and its response to local heating, and by the 
applicator form rather than the treatment duration and the irradiation 
power. 

1 What the heat diffusion limited tissue coag- 
ulation is: a phenomenon and a theoretical 
problem 

Thermal coagulation of living tissue caused by local heating due to laser light 
absorption is one of the novel thermotherapy techniques of tumor treatment 
which is currently undergoing clinical trials (see, e.g., B). Thermal coagula- 
tion is used to form a necrosis domain of desired form for the removal of the 
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Figure 1: The necrosis growth due to local thermal tissue coagulation limited 
by heat diffusion in the surrounding undamaged tissue. 



malignant tissue. So the mathematical modeling of the necrosis growth is re- 
quired, first, to find out the physical limitations and the basic characteristics 
of the treatment and, second, to optimize the therapy course. However, living 
tissue is extremely complex in structure, thereby, for the adequate theoretical 
model to be developed and for the mathematical modeling of the given process 
to be implemented reliably typical limit cases should be singled out and studied 
individually. Such a case is the subject of the present paper. 

When the tissue region affected directly by laser light is sufficiently small is 
size, the necrosis formation is mainly governed by heat diffusion. Namely, in this 
case we actually deal with the following physical process (Fig. Q). Absorption 
of laser light delivered into a small internal region of living tissue causes the 
temperature to attain such high values (about or above 70 °C) that lead prac- 
tically to immediate thermal coagulation in this region. Heat diffusion into the 
surrounding live tissue causes its further coagulation, giving rise to the growth 
of the necrosis domain. In this case heat diffusion plays a significant role in 
the necrosis growth because the necrosis size 3? exceeds the depth of laser light 
penetration into the tissue. Therefore, the temperature distribution inevitably 
has to be substantially nonuniform and for the tissue to coagulate at periph- 
eral points heat diffusion should cause the temperature to grow at these points. 
The latter property singles out the specific mode of thermal coagulation under 
discussion from other possible types of thcrmotherapy treatment and that is 
why we refer to the necrosis growth under the given conditions as to thermal 
tissue coagulation limited by heat diffusion. In particular, as will be seen below, 
the optimal implementation of this thermotherapy mode is characterized by the 
formation of necrosis domains of size 5i ~ 5-10 mm and by the treatment du- 
ration of ^course ~ 5i 2 /-D ~ 2-8 min (where D w 2 • 10~ 3 cm 2 /sec is the tissue 
temperature diffusivity) . Laser lighted is considered to be absorbed within a 
layer of thickness less than 3?. 
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Heat diffusion in the live tissue is affected substantially by blood perfusion 
causing the heat sink ]|, 0]. Thus, the characteristics of the spatial distribution 
and the dynamics of the blood perfusion rate should have a substantial effect 
on the necrosis growth limited by heat diffusion. Therefore, in modeling mode 
one has to take into account the tissue response to the temperature growth 
which can locally give rise to a tenfold increase in the blood perfusion rate Q . 
The latter effect, in particular, is responsible for a substantially nonuniform 
distribution of the blood perfusion rate and visually manifests itself in a red 
ring ( "hyperemic ring" ) appearing around the necrosis region during the ther- 
motherapy treatment. Besides, when heated the living tissue will inevitably 
exhibit spatial nonuniformities in the temperature due to the vessel discreteness 
[|| . Because of the extremely strong temperature dependence of the thermal co- 
agulation rate such temperature nonuniformities can substantially perturb the 
necrosis form and, so, this effect should be also taken into account. 

In the previous papers ||, 0, [|, [| [llj [l^, 13, basing on the free 
boundary description we have developed a model for the heat diffusion lim- 
ited thermal tissue coagulation and studied the properties of the corresponding 
necrosis growth. This model has been developed in order to take into account 
the effects mentioned above. The aim of the present paper is to survey the ob- 
tained results, to outline the general properties of the necrosis growth under the 
aforementioned conditions, and to justify the point of view on the laser induced 
heat diffusion limited tissue coagulation as a individual mode of laser therapy. 

However, before stating the free boundary model let us, first, recall the main 
features of living tissue affecting heat transfer, discuss the approach previously 
used by other authors and analyze the corresponding problems in the theoretical 
description of the local thermal coagulation. In this way we will make the key 
points of the theory we have developed more clear. 



1.1 Background: the main features of bioheat transfer 
1.1.1 Living tissue as a heterogeneous medium 

Blood flowing through vessels in living tissue forms paths of fast heat transport 
and under typical conditions it is blood flow that governs heat propagation on 
scales exceeding several millimeters (for an introduction to this problem see, 
e.g., ||, [T5J ) . Blood vessels make up a complex network being practically a 
fractal. The larger is a vessel, the faster is the blood motion in it and, so, 
the stronger is the effect of blood flux in the given vessel on heat transfer. 
The blood flux in the smallest vessels, capillaries, practically does not affect 
heat propagation [||. Thus, there should be blood vessels of a certain length 
l v that are the smallest among the vessels wherein the blood flux affects heat 
transfer remarkably. The value of i v may be different at various points (under 
nonuniform blood perfusion) and can be found from the expression |l6| (see also 
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§@): 



Here k is the thermal conductivity of the cellular tissue, c, p are the density 
and heat capacity of the tissue, j is the blood perfusion rate (the volume of 
blood going through tissue region of unit volume per unit time) and j v is its 
value averaged over scales of order of £ v (r). The factor / < 1 accounts for 
the counter-current effect. Initially it was phcnomcnologically introduced in 
the bioheat equation as a certain renormalization of the blood perfusion rate 
\ \L7\ |is| , [l9[ and its theoretical estimate will be discussed below in the present 
section. The factor L n playing a key role in the theory of bioheat transfer is 
given by the formula: 

L n = , (2) 

where l/a is the mean ratio of the individual length to radius of blood vessels 
forming peripheral systems of blood circulation. Since expression (^) estimating 
the length £ v contains the blood perfusion rate j v itself averaged over the scale 
£ v , it is implicit and has to be completed with an expression relating the value 
£ v (r), the averaged perfusion rate j„(r) and the true one j(r). At first, we write 
it in a symbolic form 

C{£ v ,j v ,j} = 0. (3) 

Nevertheless, to estimate the length £ v we may ignore the difference between the 
averaged and true blood perfusion rates in expression (|l|). Then for the typical 
values of the ratio l/a ~ 20-40 po| , the thermal conductivity n ~ 7 • 10~ 3 
W/cm-K, the heat capacity c ~ 3.5 J/g- K, and the density p ~ lg/cm 3 of the 
tissue, as well as setting the blood perfusion rate j v ~ j ~ 0.3 min -1 and the 
factor / ~ 0.5 we get from (Q) and @: 

£ v — 4 mm and L n w 3 — 4 . (4) 

The theory of heat transfer in living tissue (bioheat transfer theory) starts 
from the system of microscopic equations governing heat propagation in the 
cellular tissue and with blood inside the vessels individually. However, on one 
hand, the complex structure of the vascular network and, on the other hand, 
the lack of the detailed information about the arrangement of individual vessels 
necessitate the development of the macroscopic description of heat transfer in 
living tissue regarded as an effective continuum with certain, may be, anoma- 
lous properties. Moreover, the vascular network can vary in particular details 
from tissue to tissue or even from patient to patient for one tissue. In this case 
the development of an adequate macroscopic theory is the only way to model 
a thcrmotherapy treatment. The macroscopic governing equation for the tissue 
temperature can be obtained by averaging the microscopic equations over cer- 
tain scales that, first, are small enough so the tissue temperature does not exhibit 
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great variation on these scales. Second, these scales should be large enough for 
the heat exchange between the cellular tissue and blood could be treated within 
a continuous approximation. This heat exchange is directly controlled by the 
vessels of length £ v . So the length i v is a natural scale of the averaging in the 
theory of bioheat transfer. Besides, it turns out that the characteristic scale of 
the temperature variations is 

^ (5) 

The latter estimate demonstrates us that the value 1 /L n plays the role of small 
parameter in the bioheat transfer theory, exactly which has made it feasible to 
construct a mean field approximation regarding living tissue as a continuum 

& 

In particular, averaging the microscopic equations in this way, we have got 
the following macroscopic bioheat equation relating the tissue temperature T v 
averaged over scales of order of £ v and the averaged perfusion rate j v : 

dT v 

c P~gf = v («effVT„) - fcpj v (T v - T a ) + q v , (6) 

where we have ignored the difference in the density and heat capacity of the 
tissue and blood, respectively, c ~ c&, p ~ p bl the blood temperature T a in 
large arteries of systemic circulation is considered to be fixed, and q v is the 
heat generation rate averaged over the same scales. The factor /, accounting 
for the counter-current effect, and the ratio F = K e s/K are estimated by the 
expressions: 

f„ -L „ .5, F = F(L„)~2. (7) 
V J-'n 

It should be pointed out that the factors / and F are determined solely by 
the geometry of the vascular network, which directly results from the vascu- 
lar network being fractal in structure Q Moreover, due to the logarithmically 
weak dependence (§) these factors can be approximately treated as constants 
determined beforehand. 

It should be noted that blood flow through the large arteries can give rise 
to an anomalously fast heat transport over scales much greater than £j>. Under 
certain conditions it may appear that the effective thermal conductivity n c g 
exceeds the thermal conductivity of the cellular tissue by tenfold due to this 
effect. However, such a fast heat transport cannot be described in terms of 
the mean field theory and deserves an individual consideration. Besides, on the 
average, its role is not too essential because of the sufficiently small relative 
volume of the large arteries. Therefore dealing with local thermal coagulation 

The theoretical estimate of / ~ 0.5 has been obtained for the first time in book 1£ . 
announced in our papers [H, H] and has been later presented independently in paper 21 . 
The effect of blood perfusion on the heat propagation rate in terms of the heat conductivity 
renormalization was quantitatively studied actually for the first time in paper |22| 
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we can ignore this effect at the first approximation. We also point out that 
equation ra) is similar in its form to the phenomenological generalized bioheat 
equation p3j except for the latter contains the true perfusion rate j rather than 
the averaged one j v . This fact can be essential for substantially nonuniform 
blood perfusion jl6| as it is the case in local thermal coagulation || ^ . 

The mean field approximation leading to equation (^|) does not allow for 
the temperature nonuniformities caused by the vessel discreteness Q because 
it regards living tissue as a certain continuum. So, to take into account these 
nonuniformities we have to go beyond the scope of the mean field theory. Since 
the particular details of the vessel arrangement on scales about several millime- 
ters are practically unknown and, moreover, alter in various tissues and, may be, 
at different points of one tissue it is reasonable to regard the vessel arrangement 
and the corresponding temperature nonuniformities as random j23] . In particu- 
lar, under uniform heating, i.e. when T v is constant, these nonuniformities will 
be characterized by the mean amplitude a and the correlation length A: 

(6T(r,t)5T(r,t)) = a 2 and (8) 
(6T{r,t)ST(r',t)) < a 2 for |r-r'|>A, 

where the symbol (...) stands for averaging. As might be expected, the cor- 
relation length A is about £ v (A ~ t v ) and depends on the local value of the 
perfusion rate. By contrast, due to the fractal structure of the vascular network 
the ratio of the amplitude a to the overheating T v — T a is estimated as [[l6| : 

° < -1 (9) 

-*-v ± a 

and can be regarded as a predetermined constant. Taking into account addi- 
tional numerical factors [|l6| we get a ~ (0.1-0.2)(T — T a ) (see also ||). Since 
the vessel discreteness effect is mainly caused by the vessels of lengths about £ v 
it can be described in terms of a bioheat equation similar to (|^) which, however, 
contains the true tissue temperature T and where the averaged perfusion rate 
j v is replaced by term: 

3v->jv + 6j(T), (10) 
involving the random component Sj(r) meeting the conditions: 

(5j(r)) = 0, (11) 

(Sj(r)Sj(r')) = j 2 g° ( . (12) 



Hear the function g°(x) is such that |<7°(x)| ~ 1 for x < 1 and g°(x) <C 1 for 
x > 1. 

In this subsection we have discussed the characteristic features of living tissue 
as a heterogeneous medium. However it is also an active medium because it 
responds to temperature variations by increasing the blood vessel radius, which 
gives rise to an increase in the blood perfusion rate j. The blood perfusion rate 
can locally grow by tenfold B, so, this effect is significant and is the subject of 
the next subsection. 
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1.1.2 Living tissue as an active medium 

Living tissue tries to keep its temperature from exceeding a certain vital bound- 
ary T + fa 44-46 °C (T < T+) to prevent thermal damage (see, e.g., ||). So, 
the main increase in the blood perfusion rate j should fall on the temperature 
variations from T a to a certain value T VI ps T + and after the temperature ex- 
ceeds the value T vr the blood perfusion rate j is likely to depend weakly on 
temperature because the blood vessels exhaust their ability to expand. 

In order to analyze the temperature evolution in living tissue the dependence 
j(T) obtained experimentally under practically uniform heating of large tissue 
regions is typically used. However, whether such a dependence holds when the 
temperature distribution becomes substantially nonuniform is a question. In- 
deed, in this case due to the vessel extent the dependence j(T) can alter not only 
its particular form but also become functional. In other words, when the tissue 
is heated substantially nonuniform the blood perfusion rate j(r) at a given point 
r is determined, in general, by the whole temperature distribution {T(r')} over 
a certain neighborhood of this point (concerning a similar behavior of the tissue 
response to variations in CO2 concentration see, e.g., pQ]). Specific details of 
the mechanism by which living tissue responds to local temperature variations 
on scales about 1 cm is also a question [^0|. However, under strong heating 
localized on such scales thermal selfregulation can be effectively implemented 
through the response of microcirculatory bed to reduced O2 partial pressure or 
increased concentration of metabolism products (e.g. CO2, H + , ADP, etc.) be- 
cause higher temperatures result in a metabolism intensification with a higher 
O2 consumption. In addition, this explains a possible delay of the tissue re- 
sponse to temperature variations and enables the corresponding delay time to 
be estimated at several minutes (see, e.g., [p0|). 

In order to describe the local thermoregulation we need to specify how each 
vessel responds to the corresponding piece of information characterizing the 
state of the tissue, in particular, its temperature. So, in deriving the equa- 
tion relating the blood perfusion rate j(r,i) to the tissue temperature T(r,t) 
we make use of the following physiological data Q . The local selfregulation 
of blood perfusion in living tissue on scales about one centimeter is due to 
expansion or contraction of blood vessels making up a single microcirculatory 
bed and blood redistribution over this vascular network is mainly controlled by 
a large group of arteries differing in length significantly. The reaction of the 
microcirculatory bed is governed by receptors responding to variations in CO2 
partial pressure, H + concentration, or other metabolism products. However, 
as mentioned above, under strong heating such receptors can play the role of 
effective thermoreceptors supplying the microcirculatory bed with the needed 
information. Obviously, none of the vessels receives the whole information on 
the tissue state, so there must be a certain cooperative mechanism of informa- 
tion self-processing by which the behavior of different vessels is consistent with 
each other in the manner enabling the tissue to respond properly. 

We have shown |2?j, ^8) that such a cooperative mechanism of self- 

regulation can be implemented through the vessel response to the blood tem- 
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perature in the corresponding veins. The receptors mentioned above are located 
directly in the cellular tissue as well as embedded into the walls of vessels, in- 
cluding veins ^0|. Those embedded into the vessel walls are able to measure 
the concentration of the metabolism products directly in blood and, thus, to 
measure effectively its temperature. For small vessels (arteries and veins) the 
position of the receptors governing their behavior is not a factor. This allows 
us to make use of the proposed cooperative mechanism of self-regulation in the 
description of tissue response to local and strong heating. 

Under the adopted assumptions it turns out that for the normal tissue the 
dependence of blood perfusion rate j on the tissue temperature T can be approx- 
imately described by a local equation relating the values of j(r) and T(r) taken 
at the same point r until the temperature comes close to the vital boundary T + . 
This is due to the cooperative mechanism of selfregulation which involves the 
individual response of each vessel to the corresponding piece of information and 
the coordination of their behavior by the selfprocessing of information. The ad- 
equate selfprocessing of information is implemented through heat conservation 
during blood motion inside relatively large veins of the microcirculatory bed. 
The vascular network whose behavior is governed by this mechanism can supply 
each region of the cellular tissue with blood at such a rate that is required of its 
individual demand and different regions of the cellular tissue do not interfere 
with one other. The obtained equation for the tissue response is of the form 



where t v is the delay time of tissue response and jo is the blood perfusion rate 
provided the tissue is not affected. This result practically holds for living tissue 
containing also a small necrosis domain |l6|, [2^] because the temperature of 
blood in veins whose length exceeds its size by several times is not sensitive to 
the presence of necrosis domain. 

In the next subsection we analyze the models of local thermal coagulation 
used previously by other authors and discuss the basic theoretical problems met 
in such approaches. 

1.2 Distributed description of thermal coagulation 

The effect of heat diffusion on local thermal coagulation has been considered 
and numerically simulated by a number of authors (see, e.g. |3^, |3l], [32| and 
references therein). Leaving aside the description of laser light propagation 
in the tissue we can generalize their models for the necrosis formation to the 
following coupled equations for the temperature field T(r, t) and the field £(r, t) 
determining the fraction of undamaged tissue at a given point r and time t: 



\eee 





ja 



(13) 



dT 



V (KcffVT) - fc bPbJ {T -T a ) + q 7 



(14) 




(u,(T) . 



(15) 
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The thermal coagulation rate w(T) depends strongly on the tissue temperature 
T and for typical values of temperature attained during the treatment course 
can be approximated by the expression 



oj(T) = ujq exp 



T-T 



(16) 



Here ojq — cj(Tq), where To is a certain fixed temperature and A is a constant. 
Expression (|l^) is actually a convenient approximation of the Arrhenius de- 
pendence oj(T) ex exp{— ^} inside a not wide temperature interval governing 
the tissue coagulation. The available experimental data ]33| for the tempera- 
ture dependence of the exposure time enable us to estimate the value of A as 
A ~ 3 - 5°C (A ~ 3.6 °C for pig liver at T = 65 °C). Below in numerical 



calculations the dependence (16) will be taken in the form 

"T — 601 



u(T) = 0.2 exp 



3.6 



mm 



-l 



(17) 



where the temperature T is in degrees Celsius. 

In what follows the system of equations (0), ([H]), and expression (|l6]) will 
be referred to as the "distributed model" for local thermal coagulation. This 
model in fact is rather a phenomenological one and should be modified in order 
to allow for sophisticated effects revealing themselves in local coagulation. 

Let us consider the properties of local thermal coagulation that, in fact, 
make the mathematical description of this process stubborn and give rise to the 
inconsistency of the distributed model. During a typical course of thermothcr- 
apy treatment high temperatures of order of T max ~ 100°C are attained at 
the necrosis center, whereas at distant points T — T a w 37 °C. Under such 
conditions the damage rate u>(T) or its reciprocal value t cg (T) — l/u>(T) called 
the threshold exposure time of thermal coagulation varies extremely sharp in 
space on scales of the necrosis size !ft. Indeed, at points where the temperature 
attains, for example, the values of T = 60 °C, 65 °C, 70 °C, and 75 °C the val- 
ues of t cg are 8 min, 100, 20, and 5 sec, respectively |||. The typical duration 
^course of the thermotherapy course is about several minutes and estimated by 
the inequality, t COUISC > r [§, || , where 

t d = j^ry — 3 - 6 min, (18) 

is the characteristic time scale of the necrosis growth and (j) is the mean value of 
the blood perfusion rate attained near the necrosis boundary. In estimate ( fL8|) 
we have set / = 0.5 and used the typical values of the blood perfusion rate 
(j) ~ 0.3, 0.5, and 0.7 min -1 for stomach, intestine, and spleen, respectively J34|. 
So the layer wherein the thermal coagulation is under way at a given moment 
of time is characterized by a narrow temperature interval of thickness A located 
in the vicinity of the temperature of T cg ~ 65 °C. This layer L p( j of partially 
damaged tissue separates the necrosis region ((< 1) and practically live tissue 
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(( ~ 1 ). Thus, in order to describe the necrosis growth due to local thermal 
coagulation in the framework of the distributed model one has to consider in 
detail the temperature field inside the layer L p d, because exactly it governs the 
necrosis growth directly. 

Now we estimate the thickness <5 p d of this layer. As follows from the nu- 
merical analysis |6|, || the temperature distribution is characterized by a single 
spatial scale, so, in the necrosis domain the mean temperature gradient is about 
(VT) <~ (T max — T a )/5R. In addition, the necrosis size attaining typically the 
values of 10 — 20 mm meets the inequality [|[ |): 

SR>i T ~10mm (19) 

by virtue of (§) and (§). Then setting S pd (VT) ~ A we get the desired estimate: 



(T max Ti) 



0.24 -1mm. (20) 



So the layer L p( j of partially damaged tissue is sufficiently thin, its thickness 
is less not only than the mean necrosis size 3? but also, what is essential for 
the theory development, than the characteristic length £ v of the vessels directly 
governing the heat exchange between the cellular tissue and blood. Therefore, 
in order to describe rigorously the temperature distribution on scales about <5 p d 
we have to go beyond the framework of the mean field approximation. However, 
the distributed model is, in fact, a mean field one, because it regards living tissue 
in terms of an effective continuum. Thereby the distributed model in its own 
right cannot ensure the reliability for the mathematical modeling of local tissue 
coagulation. 

In particular, as an aspect of this problem we note that the mean field ap- 
proximation ignores the difference between the true temperature distribution in 
the tissue and the temperature field T(r, t), i.e. does not take into account the 
random component of temperature distribution caused by the vessel discrete- 
ness. This is justified if the mean amplitude a of such temperature nonuniformi- 
ties is sufficiently small. By virtue of (Ji|) and setting the temperature inside the 
layer L pd of partially damaged tissue equal to T cg « 65 °C we get the estimate 

< pd ~ A - (21) 

In other words, inside the layer L p d the amplitude of the temperature nonuni- 
formities due to the vessel discreteness turns out to be about the characteristic 
temperature drop across the whole layer L p d. So such temperature nonunifor- 
mities are bound to affect the necrosis formation. 

Another problem inherent in the distributed model is the specification of 
the main kinetic coefficients characterizing the partially damaged tissue as a 
continuous medium that effectively approximates the real cellular tissue with 
the embedded vascular network, with both of them having undergone partial 
damage. Namely, we should specify the effective thermal conductivity K e ff and 
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the factor / as functions of the blood perfusion rate j and the fraction C of 
undamaged tissue: 



Kcff = K e ff((J) , 



/ = /(C,j)- 



(22) 



However, these coefficients being the characteristics of the mean field theory can 
adequately describe heat transfer in living tissue only on scales about or greater 
than £ v . So we may write the corresponding relationships only for the live tissue 
where C ~ 1- (Clearly inside the necrosis domain K G ff = k and we may set / 
equal to any value of order unity because of j — in this region.) Inside the 
layer of partially damaged tissue any dependence of the type K e g = K ff(C) or 
/ = /(C) is nothing more than a formal phenomenological approximation of a 
real complicated phenomenon. In particular, as follows from the results to be 
obtained the particular details of the dependence / = /(C) are of no concern 
and below we will treat the value / as a constant. 

Besides, in order to complete the distributed model one should describe how 
the blood perfusion rate depends on the undamaged tissue fraction C and the 
temperature. It has been proposed [pl[ to set 



where jt{T) is the perfusion rate that would occur in tissue without damage. 
Clearly, it is also a pure phenomenological approximation. 

The only fact, that is able to justify the feasibility of applying the distributed 
model to the analysis of local thermal coagulation can be the independence of 
the necrosis growth from the particular details of heat transfer in the layer of 
partially damaged tissue. In this case, however, it would be more consistent to 
use a free boundary model which ignores the thickness of this layer, i.e. treats it 
as the boundary of the necrosis domain. The motion of such an interface must be 
governed by the boundary values of the temperature and its gradient. Particular 
details of heat transfer in the given layer may be taken into consideration by a 
certain collection of parameters. 

In the next section we will show that this independence is the case and 
derive the corresponding free boundary model. In particular, formally assuming 
equation ( |l4| ) to hold in the whole region under consideration we will reduce 
the system of equations (|l4| ) , (jl^) to a certain equivalent free boundary model 
dealing with integrated characteristics of the tissue damage rate inside the layer 
L p( j. In other words, the model to be obtained aggregates the particular details 
of the functions K c ff(C) an d j(CiT) to certain contains of order unity, which 
demonstrates the desired independence. 

2 Interface dynamics 

In order to obtain the free boundary description of local thermal coagulation it 
is suffice to consider a certain neighborhood Q of the partially damaged tissue 
layer L p( j (Fig. ||). The size of the region Q is assumed to be, on one hand, much 



(23) 
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smaller than the characteristic size 5ft of the necrosis domain and, on the other 
hand, much larger than the thickness 5 p a of the layer L p d- In other words, we 
analyze the case when the tissue in the region affected directly by laser light has 
already coagulated and the necrosis growth is caused by the heat diffusion into 
the surrounding tissue. To describe the dynamic of local thermal coagulation 
we introduce the interface T specified by the condition 

C(r,i)| rer = Co, (24) 

where £ ~ 0.5 is a certain fixed value, and keep track of how it moves. Then 
the necrosis growth is represented by the motion of the interface T at the normal 
velocity § n . 




Figure 2: The layer L pc j of partially damaged tissue and the local coordinate 
system. 



There are two ways of deriving the desired expression for the velocity $ r , 



One of them is as follows. Let us, first, assume equations (14), (15) to hold 
at all the points of the tissue. Then the derivation is to reduce the distributed 
model to certain equations dealing with the layer L p( j of partially damaged 
tissue in terms of an interface T whose velocity d n is determined by the local 
characteristics of the temperature field that practically do not vary over this 
layer and so remains unchanged at the nearest points of the live tissue and the 
necrosis region. This way has been implemented in our papers ^ using 
the singular perturbation technique in the small parameter A/(T max — T a ). In 
particular, we have obtained the following expression relating the velocity 
of the interface T, the value T cg of the temperature at this interface, and the 
boundary value of the temperature gradient, for example, on the necrosis side 
V ra T|p_ Q : 

K = »op^ ^(T cg ) , (25) 

| V„J |p_Q 

where the function uj(T) is the thermal coagulation rate given by formula ( [l6| ) 
and the constant 3q takes into account the formal dependence of the effective 
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thermal conductivity n c s{C) on the undamaged tissue fraction £ 



1 



So J Kcff(C)C 



c 



The other way of obtaining expression (|25|) is to state it basing on the general 
properties of thermal coagulation and without using the distributed model at all. 
Doing so we are not be able to get an expression for the constant and, thereby, 



should regard it as a phenomenological one. However, first, expression (26) is as 
rigorous as the former approach using the dependence K o ff(0- Second, following 
the latter way we show that the free boundary model based on expression (|2^) is 
more general than the distributed model, i.e. it is not restricted to the validity 
area of the distributed model. 

Let thermal coagulation be under way in a layer L p( j(£) at a given moment 
of time t and the mean temperature and the mean temperature gradient inside 
this layer be T cg and V ra T p d, respectively. Then the characteristic thickness <5 p d 
of the layer L p d can be estimated by expression 

S pd ~-^—. (27) 
V„ipd 

In fact, when the necrosis growth due to thermal coagulation is limited by the 
heat diffusion the temperature distribution is substantially nonuniform in space. 
Therefore, on one hand, in the necrosis region, where the temperature exceeds 
T cg by values greater than A: (T - T cg > A) and, thus, lj(T) > uj{T cg ), the 
tissue coagulation has to be completed. On the other hand, in the region of 
undamaged tissue the temperature is sufficiently low, T cg — T > A, and, so, 
oj(T) <C u(T cg ), thereby the tissue has not enough time to coagulate at such 
temperatures. Therefore tissue coagulation can be under way only in the region 
where |T — T cg | < A whence we immediately get estimate (p7|). After a lapse 
of the time interval t tin- ~ l/u(T cg ) the tissue in the layer L p d has to coagulate 
practically completely. This is equivalent to the displacement of the layer L p d 
over the distance (S p d- Thereby, if we keep watch on the points at which £ ~ 0.5 
then we will see that these points move at the velocity 

J thr Vpd-Z 

which exactly coincides with expression ( |25| ) within a factor of order unity. 

The other relation completing the free boundary approximation of the necro- 
sis growth reflects the heat flux conservation through the layer L p d of partially 
damaged tissue. Obviously, it can be represented in the form 

k V„T| r _ = k cS V„T| r+0 (28) 

relating the temperature gradients V ra T| r _ and V„T| r+0 at the interface T 
on the necrosis side and the side of live tissue, respectively. 
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Figure 3: The necrosis growth due to local thermal coagulation limited by heat 
diffusion. The physical system under consideration. 

Expressions (p5|), ( p8| ) are actually the essence of the free boundary model 
for the dynamics of local thermal coagulation. It should be pointed out that 
the stated free boundary approximation in its turn can be regarded as the 
initial point of modeling the necrosis growth due to local thermal coagulation. 
Indeed, all the information required of specifying the dependence $n(T) can 
be obtained based on the experimental data for the temperature dependence 
of the threshold exposure time at fixed temperature. The only parameter of 
this model that contains the particular information about the properties of heat 
transfer in the real layer of partially damaged tissue is the numeric factor 9 
of order unity, Qo ~ 1- Keeping the latter in mind we may regard the analysis 
presented in the given section as the substantiation of the fact that the necrosis 
growth due to local thermal coagulation limited by heat diffusion is insensitive 
to the particular details of heat transfer inside partially damaged tissue. Thus, 
it actually justifies the distributed model rather than the free boundary model. 

Now the obtained results (p5|), J2S|) enable us to state the desired free bound- 
ary model for heat diffusion limited thermal coagulation. 

3 Free boundary model 

Following the previous section the tissue is considered as comprising two regions 
(Fig. ||), the necrosis domain Q n where the blood perfusion rate is equal to zero, 

i(r,t) = for reQ n , (29) 

and the undamaged tissue Qt responding to temperature variations by increas- 
ing locally the blood perfusion rate j(r, t). Besides, we allow for that the tissue 
response can be delayed. Blood vessels can expand only to a certain extent as 
the temperature grows. So when it becomes high enough, T > T vr , the blood 
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perfusion rate j(v,t) attains a large but finite value j ma x and remains approxi- 
mately constant. Taking into account expression (|l^) we describe this behavior 
of the undamaged tissue by the equation: 

+M T )=io for reQ t . (30) 

Here t v is the delay time of the tissue response and the function $(T) is of the 
form 

s>m = / 01 + ^ ~ ^T^k for T < Tvr rsi ^ 

1 J \ a for T > T vr V ; 

where a = jo/ jmax is the ratio of jo and the maximum j max of the blood 
perfusion rate that can be attained in living tissue due to the vessel expansion 
caused by the tissue response to temperature increase, and T vr w 45-46°C is the 
temperature at which the blood vessels exhaust their ability to expand. Inside 
the necrosis domain Q n the tissue temperature obeys the heat diffusion equation 
for solids: 

cp— = k\J 2 T + q (32) 

where K is the cellular tissue conductivity and q is the rate of heat generation 
due to the laser light absorption. Inside the undamaged tissue the temperature 
is governed by the equation taking into account also the effect of the vessel 
discreteness (replacement (|To|)): 

dT 

cp ~dt = Fkv2t ~ fCbPb{jv + 6j){T - T ^ + q - ( 33 ) 

Here, as before, j v is the regular component of the blood perfusion rate averaged 
over spatial scales of order l v , the value Sj is its component due to the vessel 
discreteness treated as random, and the constants F = K e g/K > 1 and / < 1 are 
of order unity. At the interface T between the necrosis domain and the undam- 
aged tissue the temperature and the heat flux are assumed to have no jumps, 
i.e. the temperature distribution meets the following boundary conditions 

7lr+o = T\r-o = T ^ . ^V„T| r+0 = V„T| r _ . (34) 

By virtue of ( |l6| ) and ( p5|) the normal velocity of the interface V is given by the 
expression 



v n = ,„ m , exp 

|V„T| r „ 



T cg — To 



A 



(35) 



For points distant from the necrosis domain 

TL = T a . (36) 
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Finding the relationship between the averaged and true blood perfusion 
rates, j v and j, we have to take into account that the scale £ v of averaging 
in its turn depends on the local value of j v as demonstrated by expression (|I|) . 
This dependence enables us to specify functional (|^) in the form jl6| : 

j v -—^3v=3, (37) 
cp 

where A„ ~ 1/ \[T^ l is also a constant of order unity. Equation ([57]) should be 
subjected to a certain boundary condition at the interface T because it makes no 
sense to average the blood perfusion rate over the necrosis domain impermeable 
to blood. The physical layer separating the necrosis domain and the undamaged 
tissue where the local vascular network is not damaged is complex in structure 
and contains a spatial increase of the blood perfusion rate from zero to the 
value in the undamaged tissue. In order to avoid the problem of analyzing the 
blood perfusion rate in this layer we take into account the following simplifying 
circumstance. On one hand, the typical size of the necrosis domain formed 
during a thermal therapy course and the characteristic length of temperature 
variations are of the same order about 1 cm. So, particular details of spatial 
variations in the blood perfusion rate on scale much less than 1 cm are not the 
factor. On the other hand, the damaged part of the vascular network located 
inside the necrosis domain is most probable to be made up of an artery and 
vein having supplied previously this region with blood as a whole and of shorter 
vessels formed by their branching. Indeed, the mean volume of living tissue 
falling on a single small artery (or vein) of a fixed length I is about I 3 [^0|, 
p5| . Thus, the regions where the total blood perfusion is directly controlled 
by different vessel of fixed length do not cross each other substantially and the 
architectonics of microcirculatory bed can be approximately represented as the 
binary tree embedded uniformly into the cellular tissue Therefore, the 

region containing the vascular network part in which blood flow is remarkably 
disturbed because of the necrosis formation docs not exceed substantially the 
necrosis domain. The latter enables us not to make difference between the given 
layer and the interface T and to deal with a sharp jump of the blood perfusion 
rate at the necrosis interface. The desired boundary condition imposed on the 
averaged blood perfusion rate j v must obey the law of blood conservation, which 
in this case can be written as 

J drj v (r,t) = J drj{r,t). (38) 
Qt Qt 

So in order to fulfill identity (j^) we have to set the normal gradient of the 
averaged blood perfusion rate equal to zero at the interface T: 

V„j v \ r+ =0. (39) 

We note that the adopted assumption will not hold if a large vessel goes through 
the necrosis domain. However, the probability of this event is small enough and 
this case should be analyzed individually. 
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The random component Sj of the perfusion rate allowing for the vessel dis- 
creteness on scales £ v is assumed to meet the conditions: 



(Sj(r,t)) = 0, 
(8j(r,t)6j(r',t)) ~ Jv (r, t)j v (r', t)g{v, r') . 



(40) 
(41) 



Here the correlation function g(r, r') is symmetric with respect to the exchange 
of the variables r, r' , i.e. gr(r, r') = g(r',r) and has to obey the general law of 
blood conservation: 



Qt 

j v {v,t) J dv'j v (v',t)g(T,T')=0. (42) 
Qt 

For the living tissue without necrosis it is natural to assume: 



7(r,r') = <7° 







(43) 



where the function g°(x) is such that g(x) ~ 1 for x < 1, |g(a;)| <C 1 for x 3> 1, 
and 



/" da;2; 2 ff (a;) 



Q. 



(44) 



The latter reflects the blood conservation in three-dimensional space. A typical 
form of the function g°(x) is shown in Fig. 




dimensionless distance, x 



Figure 4: Typical form of the correlation function of the blood perfusion rate 
nonuniformities caused by the vessel discreteness. 
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Figure 5: The form of the necrosis region under consideration, d = 3, 2, 1 cor- 
responds to applicator of the spherical, cylindrical, and plane form, respectively 
(in the latter case ro — > oo). 

The system of equations (p9|), (|30|), (|3^), (|33| ), (|37| ) subject to the boundary 
conditions (p^)-(^6|), ( |39| ) with expressions (|40[) , (|4l|) makes up the desired free 
boundary model. 

Keeping in mind that the distributed model was widely used let us, now, 
compare the dynamics of the necrosis growth predicted with this model and 
with the stated free boundary model under different physical conditions. In this 
way, first, we will verify the basic question as to whether the value A of the 
order of 3 — 5 °C is small enough so the ratio A/ (T max — T a ) may be regarded as 
a sufficiently small parameter for the free boundary model to hold. Second, we 
will actually justify the distributed model because, as will be seen, it does be 
independent of the particular details of heat transfer inside the layer of partially 
damaged tissue. 

3.1 Model comparison 

To be defined we simulate the necrosis growth in the tissue phantom shown in 
Fig. ||. The applicator presented by the hatched circle locally heats the tissue up 
to high temperatures of the order X& ~ 100 °C, which causes immediate tissue 
coagulation in the nearest neighborhood of the applicator. This can be the case, 
for example, because of the direct heat exchange between the tissue and the 
applicator heated to such temperatures or due to irradiation of laser light and 
its absorption in an adjacent thin layer (dotted region in Fig. ||). In the latter 
case it is typically used additional internal cooling of the applicator boundary, 
so the maximum of the temperature attained just near the applicator does 
not practically depend on heat diffusion into the surrounding tissue. So keeping 
in mind the possible vaporization control over the temperature maximum we 
will treat the value 7], as a boundary temperature fixed at the interface of a 
certain radius ro. Generalizing both these situations let us confine ourselves to 
the analysis of local coagulation assuming that: 
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- at the initial time t — the necrosis region under consideration is layer of 
zero thickness whose finite radius is ro: C( r ) = and T(r) = T a ps 37 °C 
for r > ro; 

- the following necrosis growth is governed solely by heat diffusion, i.e. we 
set q(r, t) — for r > r ; 

- at the boundary r — ro the temperature is a fixed value: T(ro) = ss 
100 °C; 

- at distant points the tissue temperature is equal to T a : T — > T a as r — > oo. 

In other words, we confine ourselves to the region r > tq where the layer 
ro < r < ro+di(t) of thickness 9?(f) represents the necrosis domain whose growth 
is directly governed by the heat diffusion only. The processes in the region r < ro 
remain beyond consideration, i.e. we ignore the real dynamics of initial coagula- 
tion in the immediate vicinity of the applicator boundary. The typical duration 
of the latter process can be estimated as l/w(T(,) <C 1 sec for the applicator di- 
rectly heating the surrounding tissue and as max[l/w(Xb), (?(ro)/(cp(Tf, — T a ))] 
for the laser applicator. In the present analysis this duration is regarded as a 
small parameter. 

Keeping in mind applicators of various forms we study the necrosis growth 
in the one-, two-, and three-dimensional tissue phantoms. The tissue response 
to temperature variations is described by equation ( |30|) and expression (^l|) with 
a = 0.1 because the tissue response is sufficiently strong and the blood perfusion 
rate can locally increase by tenfold W. For the main tissue parameters we use 



the values pointed out in Sec. 1.1.1 which give us the following estimates 



To = ~ 6 min and io = \ I ~ 10 mm 

fjo V c PJ3o 

of the characteristic temporal and spatial scales, ro and £q- We compare the 
dynamics of the necrosis growth, namely, the time dependence of the thickness 
3?(t) of the necrosis layer and the coagulation temperature T cg (i) predicted by 
the developed free boundary model with that given by the distributed model. 
In the latter case we use expression ( p3| ) and specified the required dependence 
of the effective thermal conductivity K e s on the fraction £ of undamaged tissue 
by the function 

MO/* = {F - i)C + 1 • 

Besides, for the distributed mode the necrosis interface T is specified by the 
condition Cl r gr — 0-5- For the sake of simplicity in the model comparison we 
also ignore the effect of vessel discreteness and identify the true and averaged 
blood perfusion rates. 

First of all, Fig. || demonstrates the typical form of the temperature dis- 
tribution T(r), the distribution of the undamaged tissue fraction £( r ), and the 
blood perfusion rate j(r) obtained for the one-dimensional tissue phantom in 
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10 20 

depth (mm) 



Figure 6: The typical form of the spatial distribution of the tissue temperature 
(curve 1), blood perfusion rate (curve 2), and fraction of undamaged tissue 
(curve 3). (In obtaining the present curves we have set j max = 10 jo and t v rs 
1 min and considered the one-dimensional tissue phantom.) 



the frames of the distributed model. As seen in Fig. |6| the region wherein the 
undamaged tissue fraction ( varies substantially in space is sufficiently thin. So 
on spatial scales characterizing the temperature decrease such an increase of the 
value ((r) may be treated as a sharp jump. The latter actually justifies using 
the free boundary model for the given values of the tissue parameters, namely, 
the capability for assigning to an effective necrosis boundary a certain coagula- 
tion temperature T cg and certain values (on both the sides) of the temperature 
gradient. 

Fig. |?j demonstrates the fact that the distributed model leads practically to 
the same dynamics of the necrosis growth as that predicted by the free boundary 
model for one-, two- , and three-dimensional tissue phantom (figures (a, 6), (c, 
d), and (e, /), respectively). These results have been obtained for the tissue 
phantom with the strong (j max = 10jo) and delayed (t v = 2 min) response 
to temperature variations and so actually comprises the characteristic features 
of the necrosis growth in the tissue without response as well as with a strong 
immediate response. 

The present figures not only clearly justify that the particular details of ther- 
mal coagulation in a real partially damaged layer are not the factor but also 
show that the ratio A/(T max — T a ) for A w 3 — 5 °C can be treated in fact as 
a small parameter. So the two models may be regarded formally as equivalent, 
but the free boundary one does not contain self-inconsistent elements. In addi- 
tion the free boundary model can be applied to constructing a faster numerical 
algorithm of simulating the necrosis growth because in this model we need not 
consider the thin layer of partially damaged tissue. In the free boundary model 
we deal with the temperature field only in the necrosis region and the region of 
undamaged tissue where it is sufficiently smooth. 

In this analysis we have gotten another characteristic feature of local thermal 
coagulation illustrated in Fig. |^. The free boundary approximation is rigorously 
justified provided the temperature distribution inside the layer of partially dam- 
aged tissue can be regarded as quasistationary. The latter is the case when, in 
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Figure 7: Comparison of the necrosis growth predicted by the distributed model 
(DM) and the free boundary model (FBM). The thickness 9?(t) of the necrosis 
layer and the temperature Ti{t) at the necrosis interface vs time t for heat 
sources of the plane (a, b), cylindrical (c, d), and spherical (e, /) form. (In 
numerical calculations we have set t v = 2 min, j max = 10 jo- For the distributed 
model the value Tj is specified as T(r^) at the point r t at which C( r i) = 0.5.) 
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Figure 8: The time dependence of the temperature T at the point where 
C(fi) = 0.5 for different values of the parameter A. (Curves 1, 2 correspond to 
A = 1.5 °C and A = 5.0 °C, respectively. In obtaining the present curves we 
have considered the one-dimensional phantom of the tissue without response to 
temperature variation). 

particular, the time variations St,T cg of the coagulation temperature T cg are 
small enough during the necrosis growth, which in mathematical terms may be 
stated as the condition StT cg — > as A — > 0. Fig. || demonstrates that this 
condition may be fulfilled. Indeed, the smaller is the parameter A of the tissue 
damaged rate ui(T, A), to a greater extent are smothered the time variations of 
the coagulation temperature T cg except for a short initial period of the necrosis 
growth. Besides, the given feature of the local thermal coagulation justifies, at 
least at the qualitative level, the a certain simplification of the free boundary 
model || |^| that considers the coagulation temperature T cg fixed during the 
necrosis growth. 

3.2 Fixed coagulation temperature approximation 

In this approximation expression ( |35| ) is replaced by the conditions on the tem- 
perature at the necrosis boundary T. Namely, it is assumed that inside the 
undamaged tissue the temperature cannot exceed the coagulation temperature 
T cg , i.e. 

T(r,<)<T cg for r E Q t (45) 

and at the interface T the boundary value T is either equal to the coagulation 
temperature, T = T cg , or rigorously less: T < T cg . The former case takes place 
when growing near the interface T the temperature comes closely to the value 
T cg and the interface has to move in order to keep up the boundary temperature 
T inside the interval T < T cg . In the second case the interface is fixed. Both 
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these conditions can be formally described by the expression: 

dT 
~dt 



m-T cg )(. 



f) = 0, (46) 



where the boundary temperature Ti(s, t) is treated as a function of the interface 
coordinates s and the time t. 

This approximation, except for simplifying the mathematical analysis of local 
coagulation, demonstrates the fact that the necrosis growth limited by heat 
diffusion weakly depends on the particular details of thermal coagulation. All 
the characteristics of thermal coagulation are aggregated into the coagulation 
temperature T cg and the coagulation rate co(T) is only required to depend on 
temperature sufficiently strong. 

This fact is illustrated in the following way. We study the necrosis growth 
in the one-dimensional tissue phantom within the framework of the distributed 
model modified so to take into account the effect of the blood perfusion nonuni- 



formities near the necrosis domain. In other words, in equation (14) we replace 
the true perfusion rate by the averaged one, j — > j v , related by equation ([37]) 
subject to the boundary condition (p^). As before, the effect of vessel discrete- 
ness is ignored. We keep track of the points X0.2, %o.5, xo.8 specified by the 
equalities ( = 0.2, 0.5, 0.8, respectively. In this way the dynamics of the necro- 
sis growth is characterized by the time dependence of the coordinates xo.2, £0.5 
, xq.s (in mm) and the corresponding temperatures To. 2, T0.5, To. 8 (in degrees 
Celsius). Different conditions representing various possible limiting cases are 
considered. Namely, the tissue phantom without response to temperature vari- 
ations, with immediate and delayed response is analyzed. In the first case the 
blood perfusion rate remains unchanged, j = j v = jo. For tissue with imme- 
diate response (r„ = 0) the perfusion rate can attain large values directly at 
the beginning of the necrosis growth, whereas for tissue with delayed response 
(t v > 0) this increase will occur only after a lapse of time of order t„ . 

Fig. ^ shows the time dependence of the quantities xq.2, X0.5, £0.8, an d T0.2, 
Tq.5, To. s for the tissue phantom without thermoregulation (jmax/jo = 1)) with 
immediate strong response to temperature variations (jmax/jo = 5, t v = 0) and 
for the tissue where the delay in the temperature response has a pronounced 
effect on thermal coagulation (t v = 2 min). As seen in Fig. ||(61-63) the decrease 
of the temperature To. 5 becomes slow sufficiently quickly (within approximately 
1 — 2 min) and the temperature T0.5 remains inside the interval 60 — 70 °C 
during all the time of order of 5 — 10 min which is a typical duration t tr of 
thermal treatment based on thermal coagulation. The width of this interval is 
sixth as less as the typical value of the tissue overheating, T, — T a ~ 60 °C. So 
at first approximation, the temperature To. 5 may be treated as a fixed constant 
T cg . Besides, whence it follows that the coagulation temperature T cg being a 
phenomenological parameter of this approximation can be estimated from the 
expression 

*b~«- ,47) 
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Figure 9: The coordinates 2:0.2, £0.5, ^0.8 of the points at which the undamaged 
tissue fraction ( = 0.2,0.5,0.8 (al-a3) and the corresponding temperatures 
To. 2, To. 5, To. 8 (61-63) as functions of time for different value of the parameters 
a, t v . In fig. (al-a3) the thick lines labeled with FBM are the position 3? of 
the necrosis domain interface in the free boundary model with the temperature 
coagulation T cg shown in fig. (61-63). (For fig. ol, 61 - {a = 1}; for fig. a2, 62 
- {a = 0.3, t v = 0}; for fig. o3, 63 - {a = 0.3, t v = 3 min}; (A„ = 2). 
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Indeed, under such conditions thermal coagulation proceeds practically at the 
temperature T cg and the value of l/a>(T cg ) is approximately the time it takes 
for the live tissue located before in the necrosis region to be damaged. Because 
of a strong temperature dependence of uj{T) this estimate gives us the value of 
T cg to sufficient accuracy. In addition, to make the comparison of the two model 
more clear we have used in the simulation the value T cg = 60 °C found from 
expression ( pl7| ) (for t tT ~ 5 min) rather than the value of T cg approximating the 
dependence T , 5 (t) to the best degree. 

Another characteristics of thermal coagulation is illustrated in Fig. ||(al-a3). 
Comparing the time dependence xo.s(t) with the curves "FBM" (describing the 
motion 5ft(i) of the necrosis interface in the free boundary model with the fixed 
coagulation temperature) we see that there are two stages of the necrosis growth. 
The former corresponds to the time interval (0,i cg ), where t cg ~ 1/jmax ~ 
4 min. At this stage the necrosis domain grows fast enough and in the free 
boundary model the interface T reaches its limit position 5Rn m . This saturation 
of the interface displacement is due to the temperature distribution becoming 
stationary. At the latter stage (from t cg to t tr ) the necrosis domain grows 
slowly and in the free boundary model it is fixed. In other words, the proposed 
model makes the difference of these stages more pronounced. If the treatment 
is continued, the real necrosis domain will grow further and after a lapse of 
20 — 30 minutes the necrosis domain will deviate significantly in form from that 
predicted by the given model. However such a prolonged treatment is typically 
used to produce a hyperthermia effect (without visible injury) rather than to 
cause thermal coagulation directly [Q . 

It should be noted that the existence of the two stages does not obviously re- 
sult from the dependence To.^(t) because in fact the necrosis continues to grow 
slow at the second stage too. However, as follows from the present analysis, 
these stages differ from each other not only in the necrosis growth rate but also 
according to the behavior of the temperature distribution. The slow stage is 
characterized by the quasistationary temperature distribution. In other words, 
at this stage time scales on which the size 3? of the necrosis domain increases 
substantially are much larger than time scales on which temperature distribu- 
tion becomes steady-state provided the necrosis boundary is fixed. The given 
property is caused by the exponential dependence of the damaged tissue rate 
oj(T) on the temperature. At the fast stage these scales are of the same order. 

Concluding this subsection we note that the fixed boundary temperature ap- 
proximation predicts the dynamics of thermal coagulation to a sufficiently good 
accuracy. Indeed, the motion of the partially damaged tissue layer is practically 
represented in Fig. ||(al-a3) by the region bounded by the curves £0.2 (i) and 
Xo,s(t). So the free boundary model can be regarded as adequate, if the dif- 
ference \x , 5 (t) — SR(t)| does not exceed the value |xo.s(*) ~ ^0.2 (*) I remarkably. 
As seen in Fig. ^|(al-a3) this is the case except for the beginning of the growth 
and times larger than 10 min, where we have gone beyond the fixed coagulation 
temperature approximation. 

In the next section basing on the proposed model we will study the specific 
properties of the necrosis growth limited by heat diffusion. 
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4 Characteristics of the necrosis growth 



The stated above model describes several different features of living tissue affect- 
ing simultaneously the necrosis growth limited by heat diffusion. So to penetrate 
deeper in its properties it is reasonable to analyze them individually, which is 
the subject of the present section. 



4.1 Effect of the vessel discreteness 

To clarify the essence of this effect let us confine our consideration to its qual- 
itative analysis referring interested readers to specific papers [[l0| |l2|, |l3| for a 
detailed analysis. 

The vessel discreteness manifesting itself in the random component 5j(r, t) of 
the blood perfusion rate causes random perturbation <5T(r, t) in the temperature 
field. Actually 5j is the difference between the true perfusion rate and one 
averaged over scales about t v . So the latter is also the correlation length of field 



5j(r,t) (Sec. 1.1.1 ). In this subsection we mainly consider the formation of the 
necrosis domain of size 5ft ~ It ~ £ V \/L^, which enables us, at least formally, 
assume the inequality 5ft S> £ v to hold. The exception is the formation of a 
three-dimensional necrosis domain whose growth is due to a small applicator 
or the laser light absorption inside a small region Q and it will be discussed 
below. Therefore in the case under consideration the mean amplitude a and the 
correlation length A of the temperature random nonuniformities 5T(r,t) should 
be approximately the same as for such nonuniformities under uniform heating 
and be given by formula (||) and the estimate A ~ l v . The field 5T(r, t) in its 
turn affects the form of the necrosis interface T making it randomly perturbed 
as shown in Fig. |l^. It is obviously, that the correlation length ly of these 
interface perturbations is also about £ v , i.e. 

£ v ~ £ v ~ ~~^=£t ~ ^Ljft . (48) 

The remaining characteristic of the interface perturbations that should be found 
is their amplitude or, what is the same, the thickness <5r of the layer bounding 
these perturbations (Fig. [To|) . 



In the fixed coagulation temperature approximation (Sec. 3.2) the form of the 
necrosis domain perturbed by the random temperature nonuniformities 8T(r, t) 
due to the vessel discreteness is specified by the condition: 

lT v (r,t)+5T(r,t)}\ rer =T cg , (49) 

where T v (r, t) is the regular component of the temperature field (i.e. the solution 
of equations (|3^) , ( |33| ) for 5j — 0) . In a neighborhood of the necrosis boundary 
T whose thickness is much less than the mean necrosis radius 5ft the spatial 
variations of the averaged temperature T v (r, t) can be approximated by a linear 
dependence on the spatial coordinates r: 

T„(r,t) «T cg -G(r-5R) . (50) 
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Figure 10: Physical structure of the necrosis boundary. The necrosis interface 
roughness caused by the vessel discreteness. 

Here G ~ (T max — T a )/$t ~ (T cg — T a )/lT is the temperature gradient near the 
necrosis boundary Y and r = |r| provided the origin, r = 0, is placed at the 
necrosis center. Then from (E9) and ( p0[ ) we obtain GSr ~ <r whence taking 
also into account expressions (|5|), (^) we get the following estimate chain: 

J - 1 l T „ _Lsft . (51) 



n t t 

V 

Expressions dig) and (|5l| ) substantiate, at least formally, the validity of the 
following scale hierarchy: Sr < @r < 3?- So, when describing the growth of a 
necrosis domain as a whole one can ignore the effect of the vessel discreteness 
and identify the averaged tissue temperature T v (r,t) and the true one T(r, t) 
[ p"3[ . However when the attention is focused on the particular form of a necro- 
sis domain, i.e. spatial scales smaller than the mean necrosis size are under 
consideration such necrosis perturbations are substantial. 

To justify the latter statement and to measure the effect of these perturba- 
tions on the necrosis form we compare the value Sr with the thickness Sq of the 
layer of partially damaged tissue where thermal coagulation is under way. As 
results from ( p0| ) and ( pi] ) for typical values of the tissue parameters we get 

St {T cg - T a ) 

T c L n A 2 - (52) 

Therefore, these random perturbations of the necrosis interface exceed suffi- 
ciently its physical thickness, i.e. the thickness of the layer of partially damaged 
tissue. In other words, because of the temperature nonuniformities due to the 
vessel discreteness the necrosis form should be disturbed substantially, deviating 
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remarkably from that predicted by the distributed model dealing solely with the 
averaged tissue temperature T v (r, t). 

However when, for example, the characteristic size of the applicator is suffi- 
ciently small (about a few millimeters) the temperature field becomes substan- 
tially nonuniform not only due to heat dissipation caused by blood perfusion but 
also because of the geometric factor caused by the heat propagation from the 
applicator into the surrounding tissue. In three-dimensional space the temper- 
ature field T(r) induced by a small localized source even inside a tissue without 
perfusion will decrease as 1/r. Under such conditions it is this spatial decrease 
in temperature that mainly controls the necrosis growth rather than heat diffu- 
sion into the surrounding live tissue. So the thickness 3? of the necrosis domain 
is less then £t, namely, can be about 4>|§ and the perturbations of its boundary 
due to the vessel discreteness are depressed considerably Jl2| . 

Expressions (^) and (|5l| ) are actually the main results of the present sub- 
section. In the following subsections we will focus out attention on the growth 
of a necrosis domain as a whole and, thus, will ignore the vessel discreteness 
effect. 



4.2 Effects of the tissue response to heating 

In this section we analyze how the tissue response to local strong heating affects 
the necrosis growth. To single out these effects in their own right and to simplify 
the analysis we confine ourselves to studying the necrosis growth in the one- 
dimensional tissue phantom described in Sec. 3.1 (Fig. |B|) and adopt the fixed 
coagulation temperature approximation. 

First, we consider the tissue with immediate response, corresponding to 
t v = in equation (|30|). In this case the thickness 3?ii m of the necrosis domain 
attained during the fast stage of the coagulation is determined by the stationary 
temperature distribution. In particular, for the tissue phantom without thermal 
regulation (a = 1) the blood perfusion rate is constant, j — jo, and as follows 
from equations (p2|) and (p3h: 



For the tissue with thermal regulation expression (|53"|), after the replacement 
jo ^ jmax, may be also used to estimate the value !ftii m , thus 



Kum~4/ ■ (54) 

The dynamics of coagulation under such conditions is represented in Fig. [ll] for 
different values of the parameter a = jo/ j max- 

Fig. [lljz shows the size 3i of the necrosis domain as function of time. The 
higher is the tissue response, the smaller is the necrosis domain and the shorter 
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Figure 11: The size 3? of the necrosis domain vs time (a) and the temperature 
distribution (6) for different values of a. Curves 1, 2, 3 correspond to a = 1, 
0.3, 0.1, respectively (t v = 0, A„ = 2). 



is the fast stage of the necrosis growth. The duration of the fast stage is ac- 
tually the characteristic time t cg during which the size of the necrosis domain 
attains values about 5iii m . Comparing the numerical values of the corresponding 
quantities we find that the duration of the fast stage can be estimated by the 
expression 

t cg ~ j^— , (55) 

which conforms to the general properties of heat transfer in living tissue. Indeed, 
in the given model, as follows from equations ( |32| ) and (|3^), the time it takes 
for the temperature distribution to become steady-state is about /jmax and 
the establishment of this steady-state (in reality, quasistationary) temperature 
distribution is actually the essence of the fast stage. 

Fig. pd| b shows the temperature distribution for the tissue without thermoreg- 
ulation (curve 1) and for the tissue with strong immediate response (curve 2, 
imax/jo — 10). In order to compare them with each other, lengths are measured 
in units of the corresponding necrosis size. Whence it follows that, in contract 
to the time dependence 3?(t) , the tissue response to heating practically does not 
affect the form of temperature distribution. 

When the tissue response to temperature variations is sufficiently strong 
(a -C 1) the blood perfusion rate j becomes substantially nonuniform. In this 
case the averaged perfusion rate j v differs significantly from the true one j, 
which has a certain effect on the growth of small necroses. The latter feature is 
illustrated in Fig. [ij 

Fig. [T^ q compares the dynamics of the necrosis growth described by the 
given model and by the same model where, however, equation ( |37| ) is omitted 
and the replacement j v — * j is made. Fig. |l2| 6 shows the distribution of the 
tissue temperature T, the true blood perfusion rate j, and the averaged one j v 
that occur when the tissue responds in such an intensive way. In this case, as 
seen in Fig. [ilfo , the averaged blood perfusion rate can be twice as less as the 
true one. The latter has a certain effect on the necrosis growth because ignoring 
the difference between j v and j gives a more lower estimate of the necrosis size 
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Figure 12: (a) The size 3? of the necrosis domain as the function of time described 
by the stated free boundary model (curve 1) and by the same model within the 
replacement j v — > j (curve 2). (b) The distribution of the temperature T, the 
true perfusion rate j, and the averaged one j v corresponding to curve 1 at time 
t = 2 min. {a = 0.1, r v = 0, X v = 2, T vr = 50 °C) 




time (min) (depth) /(necrosis size) 



Figure 13: (a) The size 3? of the necrosis domain vs time for the tissue re- 
sponding immediately (curve 1) and with a certain delay (curve 2, t v = 2 min). 
Curve 3 represents the -Ji(i)-dependence for the tissue without thermal regu- 
lation, (b) The distribution of the temperature, true and averaged perfusion 
rates at different time moments t — 1 min (curves 1) and t = 3 min (curves 2). 
{a = 0.2, \ v = 1) 



(Fig. |15 

Now we consider how the delay in the tissue response can affect thermal 
coagulation. This effect is remarkable when the delay time t v is comparable 
with the duration of the fast stage t cg . So we may confine ourselves to values 
i~v ~ t cg . The difference in dynamics of the necrosis domain growth for the 
tissues responding immediately (t v — 0) and with a certain delay (r„ ~ 2 min) 
is illustrated in Fig. 13 and Fig. 
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As seen in Fig. |13|a when the tissue response is delayed essentially the necro- 
sis domain, firstly, grows fast enough keeping ahead of one growing in the tissue 
with the immediate response. Then the growth of the necrosis domain is sup- 
pressed and further its form remains unchanged. It should be noted that in this 
case the saturation of the necrosis domain growth is not due to the tempera- 
ture distribution becoming stationary. Its reason is the increase of the blood 
perfusion rate after a lapse of a time about r. Under such conditions there 
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Figure 14: The temperature (curve 1) and the averaged perfusion rate (curve 2) 
at the necrosis interface T vs time for the tissue with the delayed response. 
Curve 3 shows the corresponding time dependence of the necrosis size 9?(t). 
(r„ = 1.4 min, a = 0.07, A t , = 0.5) 



is enough time for the size 3? of the necrosis domain to attain values of order 
^lim 1 ^ vVio until the blood perfusion rate increases substantially. These 
values could not be attained if the tissue response were not delayed. So after 
the blood perfusion rate increases, the following growth of the necrosis domain 
becomes impossible and the temperature Tj at the interface T must go below 
T cg and the necrosis domain has to cease to grow. This behavior of the interface 
temperature is illustrated in Fig. [l4|. It should be noted that such a saturation 
of the necrosis growth for real tissues seems to be more pronounced because a 
real necrosis continues to grow at the slow stage until the blood perfusion rate 
becomes high enough. 

Fig. |l3| j> shows the distribution of the temperature and blood perfusion rate 
at different time moments for the tissue with delayed response. As before, the 
length is measured in the corresponding values of the necrosis size in order 
to compare these distributions. As time elapses the perfusion rate increases 
due to the tissue response. In contrast to this behavior of the perfusion rate, 
the form of the temperature distribution practically remains unchanged. The 
comparison of the given result and one obtained for the tissue responding to 
temperature variations immediately (Fig. pd|b) leads us to the conclusion that 
the form of the temperature distribution occurring in tissue during the necrosis 
growth depends weakly on the specific values of the tissue parameters. This 
conclusion, in particular, forms the basis for applying variational techniques to 
analysis of the thermal coagulation due to laser-tissue interaction. 

It should be noted that the conclusion concerning the universal form of the 
temperature distribution has been made analyzing the necrosis growth in the 
tissue phantoms that differ only in the properties of thermoregulation. The 
other tissue parameters (for example, the tissue thermal conductivity k and the 
initial value jo of the blood perfusion rate) take on particular values fixed in the 
present analysis. This raises the question of whether the stated conclusion will 
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Figure 15: Thickness di(t) of the necrosis layer (a) and temperature Ti{t) at the 
necrosis interface (b) vs time t for heat sources of the plane (d — 1), cylindrical 
(d = 2), and spherical (d — 3) form. Here we consider the tissue without 
response and ro = 5 mm. 



hold if we change these parameters too. However, choosing the appropriate units 
of time and length aggregating such parameters we can rewrite the governing 
equations in the dimensionless form. Thus their particular values are not the 
factor. 



4.3 Effect of the applicator geometry 

In the previous section we have actually studied the characteristics of the necro- 
sis growth when the applicator is sufficiently large and the necrosis domain can 
be locally treated as a plane layer. If the applicator is small enough so its size 
is about or even less than the thickness of the necrosis layer the heat diffusion 
into the surrounding tissue gives rise to the temperature distribution depending 
on the applicator form. The latter in turn affects the necrosis formation. This 
effect is the subject of the present section reviewing the results obtained in || . 

As in the previous section we study the necrosis growth in the tissue phantom 
described in Sec. (Fig. |^) but applying the free boundary model in its general 



form (Sec. |3|), i.e. assuming the necrosis boundary motion to be governed by 
expression (pq) because under such conditions time variations in the coagulation 
temperature become remarkable. However, keeping in mind that the latter fact 
is the main reason of the effects under consideration we will make no difference 
between the true and averaged perfusion rates in order to simplify the analysis. 

Characteristic features of the necrosis growth depending actually on the ap- 
plicator form are demonstrated in Fig. |l5|-Fig. p"7[ Namely Fig. |l5| a shows 
how the time dependence of the necrosis layer thickness 3i(t) changes for the 
applicator of plane (d = 1), cylindrical (d — 2), and spherical (d = 3) form for 
the tissue phantom with the same properties. Here is illustrated the dynamics 
of the necrosis growth in the tissue phantom without response to temperature 
variations for ro = 5 mm. Fig. p^| a shows that for one-dimensional tissue phan- 
tom the difference between the fast and the slow stages is not well distinctive 
and the necrosis growth is under way practically throughout the whole course 
of a typical thermothcrapy procedure. For the three-dimensional tissue phan- 
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Figure 16: Thickness 9?(t) of the spherical (a, b) and cylindrical (c, d) necrosis 
layers vs time t for different values of the initial necrosis radius ro = 18 mm 
(curve 1), ro = 6 mm (curve 2), and ro = 2 mm (curve 3) (a, c - tissue without 
response, b, d - tissue with strong (j ma x = 10jo) immediate response). 



torn we meet the opposite situation. In this case the necrosis growth exhibits 
the clearly evident tendency to saturation after a lapse of several minutes and 
after 12 minutes the thickness 9?| d _ 3 of the necrosis layer can attain the value 
only somewhat less than the initial radius r of the necrosis domain. At the 
same time in the one-dimensional tissue phantom the necrosis radius exceeds 
this value by several times. In the two-dimensional tissue phantom the necro- 
sis growth is characterized by the intermediate behavior. It should be noted 
that in contrast to this the time dependence of the interface temperature Ti (t) 
(Fig. |l5|6) is practically of the same form for the three cases. 

These results prompt us that the rate of the necrosis growth in the three- 
dimensional tissue phantom should depend substantially on the initial necrosis 
radius r at least for r <§C 3?| d=1 , where Sft| d=1 is the thickness that the necrosis 
layer would attain for an applicator of plane geometry. This fact is directly 
demonstrated in Fig. |l^ showing the dynamics of the necrosis growth in the 
tissue phantom without response to temperature variations (figure a) as well as 
with an immediate strong response (figure b). At a fixed moment of time the 
necrosis layers can differ in thickness by a factor of three to four as the initial 
radius ro changes from 2 mm to 18 mm. As should be expected for the two- 
dimensional tissue phantom this dependence is smoothed (Fig. |l6|) and under 
the same conditions the necrosis thickness can increase not more than two-fold. 



Fig. 17 illustrates the dependence of the necrosis growth on the initial necro- 
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Figure 17: Thickness of the necrosis layer at the fixed moment of time t — 
8 min for different values of the initial necrosis radius ro for the two- and three- 
dimensional tissue phantom, (a - tissue without response, b - tissue with a 
strong (jxaax = 10jo) immediate response) 



sis radius in the most clear form. These curves have been obtained by fixing the 
time t — 8 min and treating the point collection {9?| t _ 8 m ;„,To} for different 
values of ro as the partition points of certain continuous curves. In other words, 
this figure shows the thickness $t(t,r )\ t=s min of the necrosis layer as a func- 
tion of its initial radius ro for the fixed moment of time under various physical 
conditions. We see in Fig. [l7] that for the three-dimensional tissue phantom 
the curve 5J| d=3 (r ) practically goes into the origin of the plane {5ft, r } as for- 
mally r — > 0. When the initial radius r becomes large enough (r > 3fJ| d=1 ) 
the curve 5ft| d=3 (ro) (as it must) tends to the value 3J| d=1 (figure b). In other 
words, until r < 5ft| d=1 it is the initial radius of the necrosis domain that di- 
rectly controls the size of the necrosis layer which can be attained during a 
typical course of thermotherapy. This property of the necrosis growth is due 
to the fact that in the three-dimensional space the temperature field governed 
by heat diffusion from a local source would remain substantially nonuniform 
(T(r) oc i as r — > oo) even though we ignored the heat sink effect caused by 
blood perfusion. So in the three-dimensional tissue phantom the temperature at 
the necrosis interface T cg inevitably has to decrease substantially as the necrosis 
layer grows. So due to the strong temperature dependence of the coagulation 
rate co(T cg ) the necrosis growth will be practically suppressed after a lapse of 
time it takes for the 5ftL =1 (t) to exceed the initial necrosis radius ro. Whence 
it follows that for real applicators of the spherical form whose size does not 
exceed several millimeters the thickness of the necrosis layer will be directly 
controlled by its size and, may be, by the depth of laser light penetration into 
the tissue because in this case heat diffusion standing along leads to the sat- 
uration in the necrosis growth. In the two-dimensional space (and more so in 
the one-dimensional space) heat diffusion tends to make the temperature field 
uniform. Thus, for the two-dimensional tissue phantom the dependence of the 
necrosis growth on the initial radius should be smoothed and the thickness of 
the necrosis layer will be finite even though we formally set ro = (Fig. [l?]). In 
particular, for the case shown in Fig. [17| the change of the initial radius from 
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Figure 18: The universality of the time dependence of the coagulation (inter- 
face) temperature Tj(t). All the curves Tj(t) obtained under various conditions 
belong to the dotted region. The upper boundary of this region corresponds to 
the one-dimensional tissue phantom without temperature response, the lower 
one corresponds to the tree-dimensional tissue tissue phantom with a strong 
immediate response (j max — 10jo). 



zero to infinity causes the necrosis layer in the two-dimensional tissue phantom 
to increase in thickness by 2-3-fold. 

The dynamics of the necrosis growth depends substantially on the physical 
conditions. However, it turns out that the time dependence of the coagulation 
temperature T cg (t) (i.e. the temperature at the necrosis interface) is practically 
insensitive to particular details of the growth conditions. This fact is demon- 
strated in Fig. [ll| We have discovered that all the curves T cg (t) obtained numer- 
ically for various growth conditions are located inside the dotted region which is 
thin enough. This property of the necrosis growth is relative to the universality 



of the one-dimensional temperature distribution demonstrated in Sec. 4.2 and 
justify the conclusion that the general form of the temperature field can depend 
slowly on particular values of the main tissue parameters. Such a fact shows the 
feasibility of applying variational methods to analysis of the necrosis growth in 
more complicated cases and is undoubtedly worth special consideration. 

4.4 Self-localization of the necrosis growth in tissue with 
a tumor due to thermoregulation 



In section 4.2 we have demonstrated that the tissue response to heating affects 
substantially the necrosis growth limited by heat diffusion. It is due to the 
caused increase in the blood perfusion rate strengthening the heat dissipation in 
the surrounding undamaged tissue, which suppresses the further heat diffusion 
into it and, as the result, the necrosis growth. According to the experimental 
data p| only the normal tissue can exhibit a significant increase in the blood 
perfusion rate in the response to the temperature growth, whereas inside a 
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Figure 19: The stationary distribution of the temperature in living tissue con- 
taining necrosis domain of the size 5ft = 0.5 for two-dimensional case: (a) 
q/cp = 10 exp(— 10y), (b)- q/cp = 10exp(— 5r), r = ^/ x 1 + y 2 ,a = 0.05, 



tumor an increase in the blood perfusion rate is depressed. This leads to some 
self-localization of the temperature inside tumor and self-localization of necrosis 



growth due to thermoregulation. The latter feature is illustrated in Fig. 19, [20| 
where results of mathematical modeling of the two and three dimensional cases 
are presented. 



5 Laser induced heat diffusion limited tissue co- 
agulation as a specific therapy mode 

One of the main parameters controlling a thermotherapy treatment is the power 
of irradiation delivered into the tissue and the treatment duration. These pa- 
rameters as well as a specific therapy mode should be chosen in such a way that 
gives rise to the necrosis of desired form. Different modes, in principle, can give 
rise to a necrosis region of the same form and, in this specific mode may 

be chosen keeping in mind, for example, its optimality and stability. Various 
thermotherapy modes are singled out by those physical mechanisms that play 
the leading role and endow them with individual properties ||| . The presented 
results enable us to regard the laser induced heat diffusion limited tissue coag- 
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Figure 20: The stationary distribution of the temperature in living tissue con- 
taining necrosis domain of the size 5ft = 0.5 for three-dimensional case: (a)- 
q/cp = 12cxp(— 5r), r = \J ' z 2 + y 2 (b)-q/cp = 20 exp(—5\ / r 2 +x 2 ) 7 a = 0.05. 
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ulation as an individual therapy mode where heat diffusion into the live tissue, 
i.e. the tissue with active blood perfusion, plays the governing role. 

To justify the latter statement let us briefly summarize the aforesaid. Namely, 
we have considered thermal tissue coagulation induced by local heating due to 
laser light absorption and limited by heat diffusion into the surrounding live 
tissue. The irradiation power is assumed to be enough high for the tempera- 
ture inside the region affected directly by laser light to grow up to values about 
100 °C. Then this process is characterized by the following. 

(I) In the considered case the necrosis growth is governed by the tissue 
thermal coagulation inside a layer of thickness about 1 mm that separates the 
necrosis domain and the live tissue and involves the partially damaged tissue. 
This layer is thin enough so the standard bioheat equation cannot be applied to 
describing temperature distribution in it because on such scales this equation 
does not hold. Therefore, the approach used previously by other authors (see, 
e.g., |3^] and references therein) which is based on the continuous description 
of heat transfer and thermal coagulation inside the layer of partially damaged 
tissue (the distributed model) is inconsistent. 

(II) This inconsistency is overcome in the free boundary model of local ther- 
mal coagulation which regards the layer of partially damaged tissue as an in- 
finitely thin interface of the necrosis domain whose motion is governed by the 
boundaries values of tissue temperature and its gradient. In addition, the free 
boundary model takes into account the perturbations of the tissue temperature 
caused by the vessel discreteness. 

(III) The necrosis growth caused by heat diffusion limited thermal coagula- 
tion depends weakly on the particular details of heat transfer inside the layer 
of partially damaged tissue. This property justifies the distributed model. Nev- 
ertheless, we think that using the free boundary model is favored because in 
addition to being consistent it singles out the characteristic features governing 
the necrosis growth. 

In particular, basing on the free boundary model we have found that: 

(IV) Heat diffusion limited thermal coagulation involves two stages, fast and 
slow. It is the former stage during which the necrosis domain is mainly formed 
and its size 5i can be estimated by expression (|I|). Expression (|5^) estimates 
the duration t cg of this stage. The latter stage is characterized by a substan- 
tially slower growth of the necrosis domain, with this slowdown becoming more 
pronounced as the applicator dimensions decrease. For the applicators of the 
spherical form whose radius does not exceed several millimeters the character- 
istic size of a necrosis domain formed during a typical course of thermotherapy 
is directly determined by the applicator radius and, may be, by the depth of 
laser light penetration into the tissue. Under such conditions even heat diffusion 
itself gives rise to the practical saturation in the necrosis growth and the fast 
and slow stages are most distinctive. 

(V) The blood perfusion rate can become extremely nonuniform in the vicin- 
ity of the necrosis domain because of the strong tissue response to temperature 
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variations. In this case in modelling thermal coagulation one should take into 
account that the temperature distribution is governed by the averaged blood 
perfusion rate rather than the true one. The delay in the tissue response on 
time scales of order t cg can affect thermal coagulation essentially. In this case, 
in particular, the duration of the fast stage is controlled by the delay time r„. 
The size 3? attained by the necrosis domain under such conditions is determined 
by the initial value of the blood perfusion rate jo rather than those (j max 3> jo) 
occurring near the necrosis domain at the second stage due to the temperature 
tissue response. 

(VI) The discreteness of the vessel arrangement causes perturbations in the 
tissue temperature, which is described by introducing additional nonuniformities 
of the blood perfusion rate regarded as random. These temperature nonunifor- 
mities it turn give rise to perturbations in the necrosis boundary whose ampli- 
tude Sr and the correlation length £r are estimated as: 



Here the factor L n w ln(Z/a) w 4 (I /a ~ 40 is the characteristic ratio of the indi- 
vidual length to radius of blood vessels forming peripheral circulation systems). 
It should be noted that, first, these boundary perturbations are remarkable be- 
cause they exceed in amplitude the thickness 6$ of the layer of partially damaged 
tissue where thermal coagulation is under way. Second, the universal form of 
these relations is due to the vascular network being fractal in structure. How- 
ever, for spherical applicators of small size (about several millimeters) the effect 
of the vessel discreteness is depressed because in this case blood perfusion does 
not affect the necrosis growth substantially. 

Concerning certain additional aspects in the theoretical description of heat 
diffusion limited thermal coagulation we can state the following: 

(VII) The free boundary model can be the basis of a faster numerical algo- 
rithm of simulating the necrosis growth because it deals with solely the necrosis 
region and the undamaged tissue where the temperature distribution T(r, t) is 
a smooth field. 

(VIII) The time dependence of the temperature in the layer where thermal 
coagulation is under way (the coagulation temperature) is practically insensitive 
to the particular details of the growth conditions and to the particular values 
of the tissue parameters. For applicators of large dimensions (exceeding 1 cm) 
the form of the temperature distribution occurring in tissue during the necrosis 
growth also depends weakly on the particular values of the tissue parameters. 
The given properties form the basis for applying different variational techniques 
to analyzing heat diffusion limited thermal coagulation. 
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